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Abstract 

A general procedure for the calculation of a class of two-loop Feynman diagrams 
is described. These are two-point functions containing three massive propagators, 
raised to integer powers, in the denominator, and arbitrary polynomials of the 
loop momenta in the numerator. The ultraviolet divergent parts are calculated 
analytically, while the remaining finite parts are obtained by a one-dimensional 
numerical integration, both below and above the threshold. Integrals of this type 
occur, for example, in chiral perturbation theory at order p 6 . 
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1 Introduction 



During the last couple of years there have been growing efforts to extend the predictions 
of chiral perturbation theory, which is an effective theory of QCD at low energies, to the 
order 0(p 6 ) of the momentum expansion. Of special interest are processes which vanish at 
0(p 2 ): examples such as 77^7r°7r° |lj or 77 — >7r°77 show that the 0(p 6 ) contributions 
can be rather large but necessary to find agreement between the theoretical predictions and 
experimental data. In going beyond 0(p 4 ) one is confronted with two major difficulties. 
On the one hand, the chirally invariant structures of order 0(p 6 ) || come with a number 
of a priori unknown constants which are to be determined from experiment or estimated 
in some phenomenological model. On the other hand, at order 0(p 6 ), the diagrammatic 
expansion involves two-loop Feynman integrals. In our paper we focus on the second 
point and discuss a class of two-loop Feynman diagrams which arise naturally in the 
context of SU(3) chiral perturbation theory, namely integrals of the sunset topology for 
arbitrary internal masses, tensor numerators and (for reasons which will become clear in 
a moment) powers of denominators. The necessity to consider the arbitrary mass case 
is a characteristic distinction between SU(2) and SU(3) chiral perturbation theory. It is 
due to the presence (in the limit of isospin conservation) of three different masses (m n , 
rriK, JTijj) in the latter, as opposed to only one in the former (m n ). As a consequence, the 
calculation of the genuine two-loop diagrams becomes essentially more complex. 

Any field theory containing four-point interaction vertices with derivative couplings 
produces sunset-type integrals with numerators at the two-loop level. This is obvious 
for self-energies, but is also the case, for example, if one studies form factors for small 
momentum transfer, such as electromagnetic charge-radii of mesons. A typical contribu- 
tion is represented by the following diagram (the wiggly line is a photon, the others are 
mesons) : 




(1) 

Taylor expanding the corresponding Feynman integral with respect to the momentum 
transfer q around q = generates additional loop momenta in the numerator, increases 
the power of one denominator, and reduces the original three-point function to a two- 
point function. The resulting integrals in D = 4 - 2e dimensions [|J are of the form 



pn p p 

jr k+p,mi r k+l,m2- r l,m,3 



with integer i, j, k, I, n and propagators Pk, m = k 2 - m 2 + irj. After a tensor decomposition 
with respect to q and expressing the remaining scalar products of k, I and p in terms of 
Pk+p,mu Pk+i,m 2 i Pi,m 3 an d the Mandelstam-type variables 

si2 = ip-lf (3) 
S23 = k 2 , (4) 
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(H) can be reduced to a combination of the following basic sunset integrals 

T aAnM (ml,mlrnt-y) = «~ D f d D k d D l (5) 

fc+p,mi k+l,m,2 l,m,3 

and more simple graphs where one of the propagators has been cancelled. The latter are 
essentially products of well known one-loop vacuum integrals and thus easy to calculate, 
whereas the former are not known analytically except for special combinations of masses. 
The calculation of the quadratic and quartic charge-radii of the light pseudoscalar mesons 
in SU(3) chiral perturbation theory to the order 0(p e ), for instance, requires the general 
nonzero mass case of T a p ini \ t \ with a + (3< 7 and ri\ <6 and illustrates that rather high 
tensors and powers of denominators actually occur in practice 0. 

Sunset diagrams without numerators (a = /3 = 0) have been studied in several papers. 
In [0] they were reduced, by direct integration in momentum space ||, to two-dimensional 
integrals, which were then evaluated numerically. A one- dimensional integral representa- 
tion, derived from a Feynman parametrization, was used in M . Another one- dimensional 



integral representation was obtained from a dispersion relation in |10|, [TTj]. For the case 



of equal masses, series expansions in p 2 were considered in [[TJ, 13 1, and for the general 



mass case, various multiple series in p and the masses m\ were given in KTl |14|. In [15 



it was shown that the sunset diagram is a solution of a sixth order ordinary differential 
equation in p 2 . 

Less is known about the above integrals T a ^ )Tll>n2>n3 for positive a, (3. It is possible to 
regard them as master integrals (i.e. two-point functions with five propagators) where two 
of the propagators are massless and are raised to negative powers. In this way, expansions 
in p 2 and 1/p 2 can be calculated by the recursive procedures described in However, 



as one goes to higher orders in the expansions, the complexity of the coefficients increases 
very rapidly. Therefore, in practice, the number of terms that can be obtained is limited, 
and may or may not be sufficient, depending on the accuracy required and the value of 
p 2 considered. 

In this paper we present a hybrid method for calculating the integrals T aj p tnun2ina with 
arbitrary nonzero masses. First, we split off a Taylor polynomial in p 2 , which enables 
us to calculate the ultraviolet divergent part in an exact analytical form. Then, using a 
dispersion relation, the remaining part, which is finite as D— >4, is expressed as a one- 
dimensional integral which can be evaluated numerically to any desired precision. It is a 



generalization of the integral representation for T 0j o,i,i,i given in ||T0| , |TT |. 

Our paper is divided into three parts: in section |2| the method is described in some 
detail, in section |3] the structure of the one-dimensional integration is analyzed with 
respect to numerical instabilities, and in section £| relations among the T aj p }nun2jTl3 are 
established, which can be used to minimize the number of evaluations of these integrals 
required in the calculation of a physical process. 



2 Presentation of the method 

An important property of diagrams of the sunset topology is that their ultraviolet diver- 
gences (both the overall divergence, and those due to divergent subdiagrams) are poly- 
nomials in the external momentum p 2 . Thus, they can easily be isolated by performing 
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a Taylor expansion in this variable. Because we consider only diagrams with massive 
propagators, infrared divergences cannot occur. 

In the following discussion, the integrals (|5]) will be written as T(p 2 ) for simplicity. 
They can be continued analytically into the complex p 2 -plane, with a branch cut along 
the positive real axis from p 2 = (mi +m 2 + m 3 ) 2 to infinity. For any p 2 not on this cut, we 
may apply Cauchy's theorem, keeping D ^ 4 to regulate possible ultraviolet divergences, 

(1-T<W) = ^L/dc^jg^TK) (6) 

c 

where the operator T extracts the first r terms of the Taylor expansion in p 2 , 

r-l f^2\3 



2^ -i fa^2 



3=0 



j\ (dp 2 y 



(7) 

p 2 =0 



and C is any contour in the complex plane that encloses the points ( = and ( = p 2 and 
avoids the cut. On the cut, T(p 2 ) can be recovered from (P) by letting p 2 approach the 
real p 2 -axis from above. 

For C we choose a path consisting of two straight lines, one along each side of the cut, 
connected by an infinitesimally small circle around the branch point £= (m\ +m 2 +m^) 2 , 
and a large circle with radius R. If the number of subtractions r is sufficiently large (this 
will be specified more precisely later), the contribution from the large circle will vanish in 
the limit R^oo, yielding 

1 00 I 2\r 

(l-TM)T(p 2 ) = — / dC (( ^ )(r {n( + »0)-T((-»Q)}. (8) 

(m-i+m,2+m 3 ) 2 

At this stage, the left and right hand sides of eq. (|8j) are still dimensionally regularized. 
However, the quantity {T(( + iO) -T(£ — iO)} = 2iImT(() is finite in 4 dimensions, and 
therefore, if r is chosen large enough, the integral on the right hand side may be evaluated 
with D = 4. Thus we decompose T(p 2 ) into 

T{p 2 ) = T A (p 2 )+T N (p 2 ), (9) 

where 

T A {p 2 ) = T^T{p 2 ) (10) 
contains all the ultraviolet divergences and will be calculated analytically, and 

i 7 (p 2 )" 



= * I Kj^'^O (ID 

(mi+m2+ni3) 2 

is finite in D = 4 and will be calculated numerically. 

In order to calculate the analytic part (0), we express the derivatives of T(p 2 ) at 
p 2 = in terms of derivatives with respect to a parameter p, by which we multiply the 
momentum p: 



(p 2 Y &T{p 2 



(d P 2 y 



d 2j T((ppy 



m d P 2 J 



(12) 

p=0 
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These can be evaluated by differentiation underneath the integrals in (||). When p is 
subsequently set equal to zero, the external momentum p disappears from the propagators, 
leaving a sum of vacuum integrals of the form 

[ d D k d D l , (k-pfd-pf f (13) 

J (k 2 - mj) n i((k + If - m|)" 2 {I 2 - m|) n3 

Such integrals can be calculated recursively for arbitrary integers a', j3', JT7], [IS], |TB|. 

Let us now consider the dispersion integral ([11]) in the special case ri\ = n 2 = n 3 = 1. 
More general cases will be obtained later by differentiation with respect to the masses. 
To determine the minimal number of subtractions r required, we need some information 
about the behaviour of T(p 2 ) as \p 2 \ — > 00, which can be obtained from general theorems 
on asymptotic expansions of Feynman diagrams |I9" |. For our purposes, it is enough to 
know that, in D = 4-2e dimensions, the large p 2 expansion of Tq^i ^i has the following 
structure: 



OO OO L OO 



„ :0 (p 2 r y J .=o (p 2 ) ct " ^ ^0 (p 2 r 



(14) 

This shows that, for e close to zero, the contribution from the large circle to (§) will vanish 
in the limit R — > 00 if we choose r > a + (3 + 2. Although the calculations are simplest when 
r is chosen minimally (r = a + (3 + 2), we are free to use a larger value, which amounts 
to shifting an additional finite contribution from T N to T A . This provides us with an 
important consistency check on the method. 

An application of Cutkosky's rules [^] gives the imaginary part of Ta ^i i 1 as a two 
dimensional integral over three body phase space (the Dalitz plot): 



ImT aA x X1 {p 2 ) = J ds 12 ds 23 s" 2 s^ 3 , (15) 

P W 

where the integration region Q(p 2 ) is given by 



A 2 2 2 1 / 2 2\f 2 2\/ 2 2\ 

Am 1 m 2 m 3 + [Si 2 — m 1 — m 2 ){s 23 — m 2 — m 3 ){si 3 — m 1 — m 3 ) 
— m 3 (s 12 — m 1 — m 2 ) — m l {s 23 — m 2 — m 3 ) — m 2 {si 3 — m 1 — m 3 ) > (16) 

with 

S12 > (mi + m 2 ) 2 , S23 > (m 2 + m 3 ) 2 , s i3 > (mi + m 3 ) 2 . (17) 
In ( |161) -(|17D, s 13 = (k +p + I) 2 depends on s 12 and s 23 as follows: 

«13 = p 2 + m\ + m\ + m 2 3 - s\ 2 - s 23 . (18) 



One of the integrations in ( fl5f ) can be performed easily; e.g., at fixed s 23 , the range of 
s\2 is A - £? < S12 < A + B with 

A = 7 ^—(p 2 {ml-m 2 3 + s 23 ) - s 2 23 + s 23 (m\ + m 2 2 + m 2 3 ) - m\(m 2 2 - m 2 3 )) (19) 
^s 23 v J 

B = - — J\(s 23 ,p 2 ,m\) A(s 23 ,m|,m|) , (20) 

^•§23 



where X(x, y,z) = (x — y — z) 2 - Ayz is the Kallen function, so that 



(p-mi)'- 



/mT ai/3)lil)1 (/) = ^ 6(p 2 -(m 1+ m 2 + m 3 ) 2 ) J 



ds 23 2^ 



(m,2+m,3) 2 



a + 



3=0 



a + 1 
2j + l 



(21) 



In general, performing the second integration gives complete elliptic integrals |jrO| . How- 
ever, it is not necessary to evaluate them explicitly here. Instead, we can insert the result 
fl21 ) into the dispersion integral (|TT|) and interchange the d( and o?S23 integrations: 



;i-rW)T aAlil)1 ( P 2 



ds 23 S^ 3 1 ^/A(s 2 3,"^2,"^3) 



(m2+m 3 ) 2 



^C ,/ P i,. V / (22) 



( v / ^23+ m i) 2 



(C - p 2 ) c c 



where 



p(0 



1 j^/a + r 



p 2 =C 



(23) 



is a polynomial in £ of degree a. The integrals in fl22|) can be recognized as subtracted 
dispersion integral representations of the one-loop scalar two-point function, 



B(mlm 2 2 ;p 2 ) = -m- D/2 J d D k 



in which one of the masses is replaced with S23: 



(l-T^)B(mls 23 ;p 2 ) = J 

(v / «23+ m l) 2 



((k + p) 2 — m 2 )(k 2 — m 2 ) 



(24) 



(p 2 ) k JA(C,m?, S23 ) 



(C-P 2 )C fc 



c 



(25) 



Thus, the numerical part T^a 1 1 i(p 2 ) is given by a sum of one-dimensional integrals 
of the form 



ds 23 s 23 ^X{s 23 ,m 2 2 ,m 2 3 ) (l - T {k) ) B(m\,s 23 ]p 2 ) 



(m 2 +m 3 ) 2 



with 



k = r — a,r — a + l, ■■■ ,r 

I = (3 - a - 1, (3 - a, ... , /3 + a-l + fc-r. 

We use the following well known expressions for B(m 2 ,s 23 ;p 2 ) ||: 

1 

B(ml,s 23 ;p 2 ) = 7-/ <ixlog(m 2 x + s 2 3(l - ^) - p 2 x(l - x)) 

e Jo 

= - - 7 + 2 - logm 2 + F(m\\ s 23 ;p 2 ) , 



(26) 



(27) 
(28) 



(29) 
(30) 



where 

F(ml;s 23 ;p 2 ) = a* log (l - j-) + x 2 log (l - — ) , (31) 
xi, 2 are the zeros of the argument of the logarithm in the Feynman parameter integral, 

2^ {p 2 + ^23 - mj ± v/AbV^^)| , (32) 



and 7 is Euler's constant. The derivatives at p 2 = that enter (|26|) can easily be obtained 
from (E3), e.g.: 



2 n x 



1 , 1 1 2 . S 23 

7 + 1 - log m x H ^ 

e mf - s 23 



lo s — 

mf 



m( + S 23 

2(m? - s 23 ) 



+ 



mf 



S23J 



, s 2 3 



(33) 
(34) 



Integrals T a)i g )ni)n2jn3 where one or more of the indices n\, n 2 , n 3 are greater than one 
are related to T Q/3 1 11 by differentiation with respect to mf, m|, m|: 



1 2 2 2 2 

a,/3,ni+l,n 2 +l,n 3 +l \jn 1 ,m 2 ,m 3 ,p 



rri ( 2 2 2 2 

7^/3,1,1,1 [m 1 ,m 2 ,m 3 ;p 



(35) 



ni!n 2 !n 3 ! (<9m?)™ 1 (dml) n2 (<9m 3 ! ) n3 

Differentiating 1 1 1 in this manner leads to vacuum integrals of the same kind as in 
CP~3f) . For the numerical part, T^a ni n , it is necessary to differentiate the integrals ([26]). 
This can be done underneath the integral sign, provided one first shifts the integration 
variable s 23 — > s 23 + (m 2 + m 3 ) 2 , in order to avoid difficulties due to the square root like 
behaviour of the integrand at the lower end point. Furthermore, the mass derivatives of 
the one-loop function B are calculated most effectively using the fact that the nontrivial 
part F essentially reproduces itself: 



dB(m 2 ,s 23 ;p 2 
dm 2 



1 



A(mf,s 23 ,p 2 



(m 2 - s 23 - p 2 )F(m\\ s 23 ;p 2 ) + 2s 23 log — 

s 23 



(36) 



3 Numerical analysis 

When we are below the physical threshold, i.e. p 2 < (mi + m 2 + m 3 ) 2 , the integrands in 



(|26 ) are smooth and their numerical integration does not present any problems. However, 
some care is needed to evaluate the integrands with sufficient accuracy. The difficulties 
come from the following regions: 



s 2 3 —> oc: The expansion of B(m 2 , s 23 ;p 2 ) in l/s 23 reads (see also JR], |IU|): 



B(ml, s;p 2 ) = 7+I— logm 2 + ^ — y](jD" ( in 

€ a =l S 23 k =0 



l) 



1 m l 1 
a ka log h Oka 

S 23 



, (37) 
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where 



<t!(<t-1)! „ , 

ak ° = k\{k + ma-ma-k-l)\ **k = 0,...,*-l 
a aa = 
b „ = 

bka = a ka {if>(a) + ip(a + 1) - ip{a - k) - if) (a - k + 1)} 
for k = 1, ■ ■ ■ , a - 1 and V(n) = - 7 + t + | + ••• + 

*W = , * • (38) 
cr(cr + 1) 

By omitting all summands containing powers of p 2 less than r, it is seen that the 
leading term of (1 -T^ r ')B(p 2 ) is of order l/s^. Consequently, using the exact for- 
mulae gives rise to large cancellations, and it is much better to replace the integrand 
by its asymptotic expansion whenever S23 > S (mi + y/p 1 ) 2 , S = const. For example, 
we obtained numerically stable results with the choice S = 4 when we truncated 
the series ( ftTD at o = 50. 

S23~m^: If 7712 + 1713 <mi, the region S23~ml lies within the integration interval and 
causes problems due to factors of [m\ - S23) in the denominator (see e.g. ( |33D f). 
Nevertheless, the singularities of the different parts exactly cancel: the integrand 
is smooth at S23 ~ m\ and can be approximated by a Taylor-polynomial or rather 
a Pade-approximant. For integrals T a> p jnij i t i with a + /3<7 and n\ < 6 a diagonal 
Pade- approximation of order (20, 20) in S23 € [0.2m 2 , 5m 2 ] turned out to be a suitable 
choice. 

•S23 ~ (mi ± y/p 1 ) 2 : Differentiation of B(ml,S23]p 2 ) with respect to m\ (or S23) brings 
down a factor X(m 2 , S23, p 2 ) in the denominator, which vanishes for S23 = (mi ±y/p 2 ) 2 . 
Since B and its derivatives are smooth at S23 = (mi + y/p 1 ) 2 and (below the physical 
threshold) at S23 = (mi - y/p 1 ) 2 , a simple Taylor-approximation (up to order 10, e.g.) 
remedies the instability. 

As a rule, cancellations are more severe when higher mass derivatives are taken or more 
terms of the Taylor series need to be subtracted. 

Above the threshold, p 2 > (m\ + m2 + m^) 2 , singularities appear in the integrand at 
S23 = (Vp 1 ~ iTT'i) 2 , where the threshold of B(mf,S23',p 2 ) is crossed. Although B itself is 
still finite at this point, its derivatives are not. Therefore, noting that B(m 2 , S23',p 2 ) is an 
analytic function of S23 in the lower half of the complex plane, we deform the integration 
contour so as to steer clear of the singular point. The simplest possibility is to choose a 
straight line, 

S23 = (m 2 + m 3 ) 2 + e- 1 V, < s' < 00 , (39) 

at some angle <6 <tt/2 with respect to the real S23-axis. Along the deformed contour 
the integrand is smooth and can be integrated easilyQ The fact that the result should 
not depend on 9 can be used as a check. Notice, that such a rotation of the contour also 

-"^This approach was also used in M, where similar problems were encountered. 
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avoids the problems at S23 ~ m\ and at s 2 3 ~ (mi ± Vp 1 ) 2 mentioned above, unless the left 
end of the integration contour happens to be a critical point itself. 

We compared a number of special cases with results available in the literature. We 
computed the first 17 Taylor coefficients of 20^,1,1,1 (mf, m|, m\\ p 2 ) and verified that they 
satisfy the five-term recurrence relation derived in ||15|| . Another check involved the three 
point function (|l|), which can be calculated analytically in the special case p 2 = (p + q) 2 = 
m\ = m\ = m\ We expanded the result around q 2 = and found agreement with 
the first 3 terms as calculated by the algorithm described in this paper. We also checked 



that our algorithm reproduces the e 
which can be found in Jl3| . 



! , e 1 and e° terms of T fi t i^ i(m 2 , m? 



1 1 
m ; p 



m 



Finally, as an example, we would like to give some concrete results for 10,3,4,1,1- We 
take the masses m x = 0.28, m 2 = 1, m 3 = 1.143328474, corresponding to the masses of the 
it, K and i] 8 mesons in chiral perturbation theory (in units of mx), and two values of 
p 2 , one below (p 2 = 1) and one above (p 2 = 9) the threshold. The parameter r denotes the 



number of Taylor terms subtracted from T 0j 3 j4i1j1 and defines T ( 
fllTD . The divergent part of T A 3 4 1 1 is given by 



A 

0,3,4,1,1 



0,3,4,1,1, 



Cf. 



T, 



A,div 
0,3,4,1,1 



1 2 , A ^2 , _2 , _2^ 1 



-— [6p z + Am{ + m z 2 + m z ) + 



6e 



10m 2 



9m 2 3 



9m 2 2 



,2\2 



(v 2 ) 3 (v 2 , 

¥-4- + 12^-4- + 367p 2 + 24 7 m 2 + 67m 2 , + 67m 2 



ml 



+36p 2 log(m 2 ) + 24m 2 log(m 2 ) + 6m 2 , log^m 2 ,) + 6ml l°g( m 3 



(40) 



The finite (order e°) part of T^ 3All) and T^ 34zll are given in the following table for 
different values of r: 



r 


^0)3,4,1,1(1) 


rpAJin /-,\ 
1 0,3,4,1,1V- 1 ) 


^0)3,4,1,1(9) 


<3S,l(9) 


5 


-0.2879782058 


-43.6425974985 


3267.7085- 38384.2371 i 


-61507.6741 


6 


-0.0318842307 


-43.8986914736 


18389.8016-38384.23712 


-76629.7672 


7 


-0.0040709644 


-43.9265047399 


33170.9116- 38384.2371 i 


-91410.8772 


8 


-0.0005618259 


-43.9300138784 


49955.0126- 38384.2371 % 


-108194.9782 




^3 n 4 ,i,i(l)= -43.9305757043 


7ffi,i,i(9) = - 58239.9656 - 38384.2371 % 



The bottom row shows the total finite part TcO^i i = ^034 A 1 +3o,3Aii; which is inde- 
pendent of r, as it should be. Below the threshold, the Taylor series converges, so that 
in cases where only moderate accuracy is required (eg., 4 digits), it may be possible to 
neglect the numerical contribution T^ 34 11 . In our application 0, however, we sometimes 
needed as many as 8 digits, due to large cancellations between different T Qj/ g iniin2)n3 's. In 
any case, T^ 3 4 1 x is indispensable above the threshold, where the Taylor series diverges. 



4 Relations among the basic sunset integrals 

The method presented in the previous section enables us to calculate analytically the 
divergent part and (for given masses and external momentum) numerically the finite part 
of each T a ^^ n ,^ n , A . Nevertheless, the question arises whether the set of T at p tnitn2>nz can 
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be reduced further, both for practical and for theoretical reasons. In fact, there are three 
sources of relations, 



- permutation symmetry, if some masses are equal (e.g. T Q /3 11 1 =Tg a lilil if mi =m 3 ), 

- integration by parts identities 
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- subloop tensor decomposition, 

the third of which we would like to exploit in the sequel. 

The most nontrivial part of a fixed T Q , / g inijn2)n3 is the integral 



d D k d u l 



(l- P ) a (k 



a /p\/3 



pm pn 2 pn 3 
fc+p,mi k+l,m 2 l,m,3 



a ft pm P^i 

k+p,mi 



/Ml /Ma 

d°l pn2 • • ■ • (41) 

k+l,m 2 1,7713 



Consider for a moment just the /-integration and decompose its tensor structure in terms 
of g^ u and which is the external momentum of the /-subloop p3 |. This yields 



{k 2 



d D l 



N 



pn 2 pn 3 

k+l,77l 2 1,7713 



(42) 



for the /-integral, where N is a homogeneous polynomial of degree 2a in k 2 , k ■ /, / 2 , k ■ p 
and p 2 . After rewriting the whole integral 



d D k d u l 



(k 2 Y~ a N 



,ni pn 2 pn-i 
k+p,77ii k+l,m 2 l,rri3 



(43) 



in terms of T's again we have reduced T Qj/ g jnijn2jn3 to integrals T a > >n ' i>n ' >n ' 3 with a' = 0, 
/3' 6 {(3 - a, ■ ■ ■ , (3 + a}, n[ = m, n' 2 = n 2 , n' 3 = n 3 or w! x + n' 2 + n' 3 <ni+n 2 + n 3 . 

This procedure can be illustrated by the following graphical picture. For each triplet 
(ni, 712,113) identify the integrals T a ^ tnitn2 ^ n3 with the points of an associated two-dimen- 
sional (a, /3)-lattice and order the lattices according to n\ + n 2 + n 3 . Then the point (a, (3) 
in the (m, n 2 , ^3) -lattice can be 'projected' to the points (0, (3 -a), ■ ■ ■ , (0, (3 + a) on the 
/3-axis and points in lower lying lattices. Analogously, by changing the flow of the external 
momentum through the diagram, one can project a point on the a-axis. 




a , 

3 


1 • 


• 


• 


2 




• 


• 


1 












• 


-1 


A 




3 



On the one hand, these projections show that any T a ^ nxm<nz can be expressed in terms 
of integrals on the positive axes (a = or (3 = and min(a, /3) > 0) and trivial one-loop- 
products, so that at most the positive axes are independent. On the other hand, they 
entail further relations for fixed (ni,n2,n 3 ): 



10 



- First, there are two possibilities to project the points (a, a) on the diagonal onto 
a positive axis, which sets up a relation between the point (2a, 0) and the points 
(a', 0) with < a' < 2a and (0, (3) with < (3 < 2a. In other words, the even points on 
one of the positive axes, e.g. the positive a-axis, are eliminated from the minimal 
set of T's. 

- Second, equating the two projections of the point (a,a + l) in the line below the 
diagonal (cf. the above figures) yields a relation between the point {2a + 1,0) and 
the points (a', 0) with < a' < 2a, (0, 0) with </? < 2a + 1 and (-1, 0). The integral 
(—1,0), which is of a different topology, can be expressed in terms of (1,0), (0,0) 
and (0, 1) by equating the two projections of (0, 1). In this way the odd points of 
the positive a-axis with the exception of (1, 0) are eliminated from the minimal set 
of T's. 

The point (1,0) can be eliminated as well by equating the two projections of (2, 1) and 
using the above results for (3, 0) and (2, 0). But, unlike the relations for the points (a, 0) 
with a > 1, the relation found in this way contains a denominator (m 2 -p 2 ){m 2 , - m^) and 
is therefore not applicable in all mass cases. For this reason, we prefer to keep T lj0jniin2in3 
in our set of basic integrals, along with the integrals To i( g iniin2in3 with f3> 0. 



5 Conclusion 

In this paper we presented a method to calculate the dimensionally regularized sunset 
integrals T Qj/ g jnijn2jn3 (m^, m|, m|; p 2 ), see (^), with arbitrary nonzero masses and arbitrary 
momentum p 2 , both below and above the threshold. The simple topology of such dia- 
grams - they can only be cut in one way - allows their divergences to be extracted by 
subtracting a finite number of terms of their Taylor expansion in p 2 around p 2 = 0. The 
subtracted terms are vacuum diagrams that are known analytically. The remainder, which 
is finite in D — 4 dimensions, is transformed into a well-behaved one- dimensional integral 
that can be evaluated numerically. Furthermore, we discussed a set of linear relations 
between different T a ,/?,rai,n2,n 3 's and suggested a reduction mechanism to a basic set of 
these integrals, namely T lj0;num>ri3 and T 0i/3iniin2jn3 with /5>0. 

The algorithm was developed in order to meet the needs of a calculation in SU(3) chiral 
perturbation theory, but it can be applied in any field theory describing massive particles 
with four-point interaction vertices. The basic idea of the method can be extended to any 
two-loop diagram containing a subgraph which itself corresponds to a two-point function 
[PTJj ; in particular, the full three-point graph ([!]) can be handled in this way for arbitrary 
momentum transfer q. For chiral perturbation theory the techniques presented in this 
paper are a step towards taking into account the full mass dependence of the nontrivial 
two-loop diagrams which occur at the order 0(p 6 ) of the momentum expansion. 

We thank D.J. Broadhurst for some valuable suggestions, J. Gasser for drawing our 



attention to ref. ||21|| , and K. Schilcher for helpful discussions. J.B.T. was supported 
by the Graduiertenkolleg "Elementarteilchenphysik bei mittleren und hohen Energien" in 
Mainz and P.P. by the "Studienstiftung des deutschen Volkes" . 
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